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Abstract 

The production of two bb pairs is a prominent background for Higgs and New Physics searches 
in various extensions of the Standard Model. We present here the next-to-leading order QCD 
corrections to the quark induced subprocess using the GOLEM approach for the virtual corrections. 
We show that our result considerably improves the prediction and conclude that the inclusion of 
next-to-leading order effects is indispensable for reliable studies of bbbb observables in hadronic 
collisions. 
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1. Introduction 

With the LHC data taking on the horizon, the quest for precise predictions for multi-parton 
processes becomes more and more pressing. The confirmation of the Higgs sector of the Standard 
Model (SM) and the exploration of the TeV scale with respect to New Physics will ultimately 
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depend on a very detailed understanding of signal and background processes [H, 0) S 01 • Especially 
in hadronic collisions, scale uncertainties lead to a strong logarithmic dependence if only leading 
order (LO) perturbation theory is applied. Whereas the inclusion of next-to-leading order (NLO) 
QCD effects is relatively straightforward for simple kinematic situations already for four parton final 
states the number of results is very limited. The complexity of next-to-leading order computations 
with many external partons are many-fold. It is not only necessary to provide a numerically stable 
implementation of the ultra-violet (UV) renormalised one-loop amplitude, one also has to combine 
the different parts of the computation such that infrared (IR) divergences are compensated between 
virtual and real emission corrections. After IR cancellations the integrands are finite and a stable 
integration over phase space can be achieved. In all respects a lot of progress has happened in 
the last few years which point towards automated approaches concerning the one-loop amplitudes 
„„„„„ ^[^[^ 2i, Hal which 
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can be combined with publicly available Monte-Carlo tools 24> 25, 2^, 23,125 

Recent achievements are the evaluation of the order as corrections to pp — > bbti 
pp — > Wjjj [3, in, m, nil . In both cases the inclusion of higher order effects leads to the expected 
improvement in predictability of corresponding production rates which is of importance for Higgs 
and New Physics searches. 

In this paper we present a result of similar complexity which is relevant in the context of Higgs 
searches in the Minimal Supersymmetric Standard Model (MSSM) and two Higgs doublet extensions 
of the SM, where five scalar Higgs bosons are present. To scan and explore the parameter space in 
these models Higgs pair production is a relevant handle. As light Higgs bosons prefer to decay into 
6-quarks in large parts of parameter space, experimental studies are so far strongly affected by the 
uncertainty of the Standard Model backgrounds, especially pp hhbb, see for example Ref. (ssj . 
Although mainly motivated by supersymmetry, the four-6 final state also allows the study of other 
interesting beyond Standard Model (BSM) scenarios such as hidden valley models, where decays of 
hadrons of an additional confining gauge group can produce high multiplicities of bb pairs 0, Is^ ■ 
This is the reason why the pp — > bbbb process was added to the experimenters wish-list of relevant 
next-to-leading order computations Q . 

For the evaluation of the process pp — > bbbb two partonic initial states contribute at tree level: 
gg bbbb and qq bbbb. In this article we present the next-to-leading order QCD results for the 
^ark induced process. Our calculation is based on the GOLEM approach for one-loop amplitudes 
43 [which has already been applied in amplitude computations of various complexity 
47|. The IR divergences are treated using the Catani-Seymour dipole subtraction 
method [49| implemented in MadDipole [g^]. For the evaluation of tree amplitude contributions 
the matrix element generators Madgraph/MadEvent [iS, llHl and Whizard (27I are used. In the 
latter the relevant dipole subtraction terms have been included, too. 

In the next section we will briefly sketch how our result is obtained and in section 3 numerical 
results are presented. Section 4 resumes the paper. 




2. Calculation 

At leading order two topologies contribute to the amplitude qq bbbb, see Fig. ([T]). As the 
parton distribution functions for the 6-quarks are very small we only consider the cases q — u,d, s, c. 
We treat the 6-quarks as massless which is a very good approximation for LHC kinematics and cuts 
as long as the &-quarks can be detected and are sufficiently separated in phase space. Effects of the 
heavy top quark are neglected altogether as they are numerically not important. 
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Figure 1: Six-quark topologies which contribute to the leading order amplitude qq — > bbbb. 

2.1. Virtual corrections 

The evaluation of the virtual corrections is the technically most challenging part of the calcula- 
tion. Our approach is based on the Feynman diagrammatic representation of the amplitude. The 
most complicated topologies are of pentagon and hexagon type as shown in Fig. We evalu- 




Figure 2: Pentagon and hexagon one-loop topologies which contribute to the virtual corrections to order ©(af ). 

ated all Feynman diagrams in two independent ways using the 't Hooft-Veltman scheme. What 
is common to both methods is the projection onto different helicity and colour structures. One 
strategy (1) is to convert the Feynman diagrams to a form factor representation and to evaluate 
these and the resulting objects numerically. In this case we use Feynman diagram representations 
generated with QGRAF jSOj ■ FOR M |5l| is used to map the helicity amplitudes to a tensor form factor 
representation as defined in (8l.ll6l.l53|. The latter is exported to FortrELn95 code such that it can 
be linked to the form factor library golem95 For the code generation we have developed a 

dedicated optimisation tool which can also be used for other purposes [sl]. The other strategy (2) 
is to perform a symbolic tensor reduction down to genuine scalar integrals. Here we use FeynArts 
and FeynCalc [52] to generate the amplitudes. The tensor integrals are reduced symbolically to 
scalar integrals also using the formalism described in Q. As a scalar integral basis we choose two- 
and three-point functions in Z? = 4 — 2e-dimensions and four-point functions in 6 dimensions. In 
this way the extraction of infrared and ultraviolet poles is very transparent. In both strategies the 
reduction of the five-point tensor integrals is done without introducing inverse Gram-determinants. 

Without UV counterterms and IR subtraction terms one has non- vanishing coefficients for double 
and single poles in 1/e. Adding in the UV counterterm and the Catani-Seymour IR insertion 
operator, J(e) [i^ . results in the numerical cancellation of the pole part. This by itself is a stringent 
check on the algebra. 

All Feynman diagrams are thus evaluated in two completely independent manners. Together 
with the cancellation of UV/IR poles and the test of symmetry properties of helicity amplitudes 
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the correctness of the virtual ampHtudes is weU checked. 

2.2. Real emission corrections 

Apart from the one-loop corrections we have to consider processes with an extra parton in the 
final state. In this paper we restrict ourselves to qq initial states. Further we are only interested 
in observables with four well separated b-jets where the 6-quarks can be tagged. The relevant 
amplitude is thus qq — > bbbbg. If the extra gluon is soft and/or coUinear infra-red divergences 
are treated by using the dipole subtraction method of Catani-Seymour _49]. In our case 30 dipole 
subtraction terms are needed. In processes with so many dipoles the efficiency of the integration 
over phase space can be improved considerably if the dipole terms are only subtracted in restricted 
regions. In the approach of [55] a slicing parameter a controls the size of phase space where the 
subtraction terms are included. This method is adopted here. Eventually the results must not 
depend on the slicing parameter a. We have carefully checked that this is indeed the case. 

For the evaluation of all tree like contributions of our process, i.e. the LO contribution and 
the real emission corrections including dipole subtraction, we use adapted versions of the tree level 
event generators Madgraph/MadEvent [24, 25] and Whizard [27l [ist. For the IR subtractions we 
use MadDipole (2(]| supplemented by integrated dipole contributions. Again we have two fully 
independent implementations which should guarantee the correctness of the result. 

3. Results 

In this section we provide our prediction for the quark induced four 6-quark production rate 
at next-to-leading order QCD for the LHC. We use ^/s = 14 TeV as center-of-mass energy for all 
results. To define a bbhh event we use the following strategy. 

1. ) First we apply a jet algorithm to identify the phase space region where the 6-quarks can be 

separated from the extra gluon jet by applying the kr jet algorithm as explained in [H^. If a 
gluon is merged with a &-quark we obtain an effective 6-quark with a momentum pj, = pg +pb. 

2. ) We apply separation cuts on the four resolved 6-quarks which are defined by the transverse 

momenta, rapidities and separation parameters as 

PT{bj) > 30 GcV 

1^(^)1 < 2.5 
AR{b„bj) > 0.8 . (1) 

The separation cut is defined in terms of the azimuthal and rapidity differences of the (effectively) 
observed &-quarks, i.e. AR{bi,bj) — y^{4>i — 4>jY + {rji — rjjY. Note that in the case of an unre- 
solved gluon/ & system the effective momentum p^ is used to test the cut conditions. We use the 
CTEQ6M parton distribution functions (13] with 2-loop running for as for both the LO and NLO 
order cross section evaluations. 

In Fig. [3] we show the renormalisation scale (/^_r) dependence of our LO and NLO results. 
We vary the renormalisation scale around a standard scale which we define as the sum of 

the squared transverse 6-quark momenta, /io — \J^j Pri^j)- '^^^ factorisation scale is fixed to 
/ii? = 100 GeV. As not all initial states are included yet we postpone the study of the ^xp dependence 
to a forthcoming paper. As expected the LO result shows a strong monotonous dependence on 
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Figure 3: The dependence of the LO and NLO prediction of pp(qq) — > bbbb + X ai the LHC = 14 TeV) on the 
renormaUsation scale = x^o with /^o = PtC'j )- factorisation scale is fixed to = 100 GeV. 



X ~ fJ-R./ IJ-0- Inclusion of the NLO corrections leads to a cancellation of the leading logarithmic 
dependence such that the NLO result has now a plateau for a value around x ~ 1/2, where the 
prediction shows the desired stability under scale variation. 

As a reference we provide cross section values for four scale choices in Tab. [T] 





<JLO [fb] 


(JNLO [fb] 


2fio 


197.1 


342.7 




286.8 


413.8 


Mo/2 


434.0 


463.0 


A'o/4 


692.4 


434.0 



Table 1: Total LO and NLO cross sections for different renormaUsation scale choices. As a standard scale we use 
HO = ^/y^^- p^(b, ). The factorisation scale is set to = 100 GeV. The precise definition of the four-6 observable 
is given in the text. The error on all these numbers is estimated to be well below 1%. 

Finally, we compare the LO and NLO prediction of the invariant mass distribution of the leading 
&-pair in Fig 21 This invariant mass distribution is defined by the leading and sub-leading b (b) 
quark, i.e. we choose the two quarks with the largest pr and evaluate the invariant mass of this 
system. We define an uncertainty band by varying the renormaUsation scale around the scale 
between /io/4 and 2/io, where again fio = Pri^j)- This asymmetric scale choice is motivated 

by the asymmetric shape of the NLO dependence, see Fig. |3l For the scale iiu — /io/2 the LO 
and NLO distributions look very similar which is indicated by the full (NLO) and dashed (LO) 
line in Fig. 2] This is consistent with Fig. [3] as for this value the inclusive cross sections almost 
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Figure 4: Invariant mass {mi,},) distribution of the two leading 6-quarks (see text). The LO/NLO bands are obtained 
by varying the renormalisation scale fin between /io/4 and 2^o with fio - 
shows the NLO (LO) prediction for the value fin = fio/2. 



Y^.p^ibj)- The full (dashed) line 



coincide. The importance of the higher order terms is evident from the uncertainty band which is 
considerably reduced after inclusion of the NLO QCD corrections. 



4. Conclusion 

In this letter we have presented the next-to-leading order QCD prediction for the quark induced 
four-6 cross section at the LHC. It is the first part of a complete NLO QCD evaluation of pp — > 
bbbb + X. This final state is an important background for BSM Higgs searches at the LHC as it is 
a prominent final state for Higgs pair production in large regions of parameter space. 

From the technical point of view we have provided another example for the Feynman diagram- 
matic approach being well capable to compute one-loop processes of such a complexity. We have 
used the GOLEM method implemented in a highly automated framework such that further results 
can be obtained with the same set-up. Two independent evaluation strategies guarantee the va- 
lidity of the result. For the tree-level amplitudes we use Madgraph/MadEvent in combination with 
MadDipole and Whizard. For the latter the relevant dipole terms have been added. We observe 
that the evaluation time for any observable is dominated by the real emission part of the calcula- 
tion. Future improvements can be expected by phase-space parametrisations which are adapted to 
a given IR subtraction method. 

We have shown for the quark induced case that the inclusion of the order as corrections leads to 
a significant improvement of predictions for four 6-observables at the LHC as the renormalisation 
scale dependence is considerably reduced. The same can be expected for the the gluon induced case 
which is under construction using the same methods. 
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